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The matrix product state (MPS) is utilized to study the ground-state properties and quantum phase transitions 
(QPTs) of the one-dimensional extended quantum compass model (EQCM). The MPS wavefunctions are ar- 
gued to be very efficient descriptions of the ground states, and are numerically determined by imaginary time 
projections. The ground-state energy, correlations, quantum entanglement and its spectrum, local and nonlocal 
order parameters, etc., are calculated and studied in details. It is revealed that the von Neumann entanglement 
entropy, as well as the nearest neighbor correlation functions, can be used to detect the second-order QPTs, but 
not the first-order ones, while fidelity detections can recognize both. The entanglement spectrum is extracted 
from the MPS wavefunction, and found to be doubly degenerate in disordered phases, where nonzero string 
order parameters exist. Moreover, with linearized tensor renormalization group method, the specific heat curves 
are evaluated and their low temperature behaviors are investigated. Compared with the exact solutions, our 
results verify that these MPS-based numerical methods are very accurate and powerful, and can be employed to 
investigate other EQCMs which do not permit exact solutions at present. 
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I. Introduction 

During the past several decades, the role of the orbital 
degrees of freedom in determining the magnetic and trans- 
port properties of transitional-metal oxides (TMOs) has been 
widely recognizedi 1 - - — The complex intrinsic interplay in 
TMOs induces their extremely rich phase diagrams and var- 
ious fascinating physical phenomena. In order to mimic 
these orbital states with a two-fold degeneracy, the quan- 
tum compass model (QCM) was firstly introduced by Kugel 
and Khomskii. 7 The orbital degrees of freedom are repre- 
sented by pseudospin-1/2 operators, and the competition be- 
tween orbital orderings in different directions is simulated 
by anisotropic couplings between these pseudospins. Par- 
ticularly, two-dimensional (2D) QCM has attracted consid- 
erable attention due to its interdisciplinary character. Be- 
sides the ability to describe t^g systems, it was also proposed 
that the compass model can describe the physics of protected 
qubits, 8,9 and hence it may have potential application in the 
quantum information techniques. The strong quantum frus- 
tration makes it difficult to solve the system analytically, and 
consequently leads to large degeneracy in the energy spec- 
trum, which sets obstacles for numerical simulations,— It is 
generally implied that there exist a symmetry-broken ground 
state and a first-order quantum phase transition (QPT) at the 
self-dual point. — ~— 

On the other hand, the one-dimensional (ID) QCM has 
also triggered extensive studies^— In Refllol by mapping 
to the quantum Ising model, Brzezicki et al. obtained an 
exact solution of the ID extended QCM (EQCM), reveal- 
ing that it exhibits a first-order transition between two dis- 
ordered phases. Subsequently, Wen-Long You and Guang- 
Shan Tian adopted the reflection positivity technique in the 
standard pseudospin representation to rigorously determine 



the ground-state degeneracy. And, a first-order phase tran- 
sition was also confirmed. Following the approach in RefJIfJ 
Eriksson and Johannesson 22 studied the QPTs in a ID EQCM 
with more tunable parameters. They suggested that the re- 
ported first-order phase transition in fact occurs at a multi- 
critical point where a line of the first-order transition meets 
with a line of the second-order transition. Generally speak- 
ing, a first-order QPT is often associated with energy level 
crossing in the ground state, and hence the entanglement mea- 
sures, such as concurrence and entanglement entropy, would 
behave discontinuouslyi 2 ^ 2 ^ However, in Ref. 22, the authors 
claimed that they encountered an "accidental" exception. The 
concurrence and block entanglement can accurately signal the 
second-order transitions, but not the first-order ones. In other 
words, the entanglement measures do not show any discon- 
tinuities or singularities across the first-order quantum criti- 
cal points (QCPs) in the ID EQCM. Nevertheless, a converse 
point of view that both concurrence and quantum discord can 
reliably detect the first-order QPTs of this model was pro- 
posed very recently. 26 In a sense the ID EQCM can be exactly 
solved by taking Jordan- Wigner transformation, nevertheless, 
it is still not easy to analytically calculate the spin correlations 
for arbitrary sites and the excited states. 

In this paper, we investigate numerically the ground-state 
properties and QPTs of ID EQCM with matrix product state 
(MPS) variational wavefunction and the related algorithms. 
We would like to point out that MPS is a very useful and 
highly efficient real-space description of the ground states, 
and it provides a novel way to study the QPTs in EQCM. To 
be specific, firstly, some exact MPS ground states for EQCM 
Hamiltonian in some limiting cases can be obtained, and for 
off-limiting generic parameters, infinite time-evolving block 
decimation (iTEBD) algorithm 27 is adopted to determine the 
variational MPS ground state. Very accurate results can be 
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achieved in gapped regions (up to 8~9 digits compared with 
exact solution, see Fig. [3] below) with a small number of re- 
served states. In addition, the iTEBD method can also be em- 
ployed to take adiabatic continuation calculations, which ap- 
parently reveals the energy level crossing around first-order 
QPTs. Secondly, given the real space wave-function in MPS 
form, the interesting quantities including the ground-state en- 
ergy, energy spectrum, correlation functions, entanglement 
entropy, fidelity per site, as well as local and nonlocal or- 
der parameters, etc, can be conveniently evaluated. Some of 
them are not easy to be obtained by other methods. Thirdly, 
the MPS-based algorithms can be applied to other extended 
models, and hence provide us powerful tools to explore other 
EQCMs without exact solutions. 

Through the numerical calculations with MPS, we verify 
the phase diagram of the ID EQCM (see Eq. (Q~|i), and it is 
uncovered that both the first- and second-order QPTs can be 
detected by the fidelity, while the entanglement measures can 
only capture the later ones. Furthermore, we discover that the 
entanglement spectra in disordered phases of EQCM happens 
to be doubly degenerate, and correspondingly there exist two 
non-local string order parameters, which reveals the hidden 
Z-i x Zi symmetry breaking. 

This paper is organized as follows. In section II, the Hamil- 
tonian of the ID EQCM is introduced, along with the MPS 
description and related perturbation analysis. Besides, the en- 
tanglement and fidelity measures in the framework of MPS 
are concerned and discussed. In sections III, we provide our 
main numerical results, which include the ground-state en- 
ergy, entanglement entropy, fidelity and string order param- 
eters in different regions of the phase diagram. Afterwards in 
section IV, with finite temperature algorithm, i.e., linearized 
tensor renormalization group (LTRG), the specific heat curves 
of ID EQCM are calculated and analyzed. Finally, some pos- 
sible extensions of present work, as well as a summary, are 
presented in section V. 



II. Model Hamiltonian, Matrix Product State, Entanglement 
Entropy and Fidelity 
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FIG. 1: (Color online) Schematic phase diagram of the ID EQCM, 
the four different phases are marked as regions I, II, III, and IV. The 
dashed lines denote three typical paths ((a) J2 = 2(1 — Ji); (a) 
J2 = 1 — Ji ; (c) J2 = 0.8(1 — Ji )) that will be used in the following 
discussions. L\ = 1 is set as the energy scale. 



A multicritical point (Jj = 0, J2/L1 = 1) locates where the 
lines of the first-order and second-order QPTs meet. 

The Hamiltonian ([T|i commutes with the parity operators 
Pi = cTzi-i&^i' an d tnus me parity of every odd bond is con- 
served. In such circumstance, the Hilbert space can then be 
decomposed into subspace V(pi), where pi is the eigenvalue 
of Pi and introduced to label the relative pseudospin direc- 
tion on odd bonds, that is, pi = (pi = 1) when the two 
pseudospins are parallel (antiparallel). It is disclosed that the 
ground state lies in space {pi = P2 = ■■■ = PN' = 0} for 
Ji < 0, and in {pi = p 2 = ... = pn' = 1} for J\ > 0, 21 
that is, for ferromagnetic coupling J\, two spins on odd bonds 
can only be parallel (one of such spin configurations is shown 
in Fig. 12(a)), while for antiferromagnetic (AF) coupling J 1; 
the spins on odd bonds must be antiparallel in the ground state 
(see Fig. |2](b)). 



A. Quantum Compass Model 
The ID EQCM is given by 

JV' 
i=l 

where periodic boundary condition is assumed, and N = 2N' 
is the total number of sites. The a\ ,z are Pauli matrices on the 
zth site, Ji, J2 on odd bonds, along with L\ on even bonds, are 
exchange couplings. For the following calculations in sections 
III and IV, the coupling constant L\ = 1 in Hamiltonian Eq. 
flj is set as energy scale. 

The ground-state phase diagram of ID EQCM (see Fig. [U 
is sketched by previous studies^ As is shown in Fig. [T] the 
system undergoes a first-order QPT identified by critical line 
J 1 = and a second-order QPT with critical line J%jL\ = 1. 



B. Matrix Product States and Perturbation Analyses 

In principle, any quantum state of many-body system can 
be expressed by MPS form through taking successive Schmidt 
decompositions site by sitei^22 However, not all of these 
MPS expressions are efficient and can be utilized for sim- 
ulations. For ID cases, thanks to the entanglement entropy 
area law, 30 quantum spin chain with only finite range interac- 
tions and possesses a gapped spectrum can be efficiently sim- 
ulated with the MPS-based algorithms. MPS is closely related 
with the density matrix renormalization group method^ and 
it well satisfies the entanglement area law in ID. When the the 
chain is divided into two blocks by cutting an bond, the renor- 
malized left (right) bases are the eigenvectors of the reduced 
density matrix subsystem to the left(right) of the broken bond. 
As long as the entanglement entropy between two subsystems 
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FIG. 2: (Color online) (a) and (b) show two typical spin configura- 
tions, the spins on odd bonds are parallel (for Ji < 0) or anti-parallel 
(for Ji > 0). The dash ovals denote the spins on even bonds, (c) and 
(d) show the one- (r) and two-period (T a , Tt) MPS wavefunctions, 
A is a diagonal matrix on each bond. 
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FIG. 3: (Color online) The calculation errors of energy per site, eMPS 
means MPS numerical result and e cx is the exact solution. For non- 
critical regions, even with smallest bond dimensions (x = 2 or 4), the 
MPS calculation provides very accurate results. The systems at or 
near the critical line ( J2 / Ji = 1.0 and 1.05 shown above) are harder 
to tackle, however, a few more states (say, x — 20 or 30) is adequate 
in practical calculations. 



is bounded by the area law, the classical simulation of the ID 
quantum system can be performed efficiently. 

For the present ID EQCM, as shown in Fig. |2](c), a tensor 
pni2i,m 2 , + i j s use( j t0 address the two spins on even bonds, 
where m means local spin physical index, and the wavefunc- 



tion in uniform (period-one) MPS form can be written as 



N' 



I* 



MPS, 



= Ti(JJr Tnai ' m2<+1 A i )|..., (TO 2l ,m 2i+ i), 



(2) 



in which, A means a x x X diagonal matrix, x 1S also called 
the bond dimension, and Tr is the trace of matrix product. 

To explain this point more explicitly, we take the limiting 
cases with L\ = into account, where the spins on even 
bonds are unentangled, and the exact MPS ground states are 
thus obtainable. For J\ < (and J 2 > 0), the limit parameter 
point belongs to region I of the phase diagram Fig. [T] The 
ground state of the local two-site Hamiltonian on odd bond is 



1^(2* -1,2*)) 



1 



V2 



I t2i-lt2i — \-2i-l-i-2i) , (3) 



with bond energy J\ — J 2 , and it requires that the spin orien- 
tations on odd bonds must be parallel. Thus MPS Eq. (|2]i with 
the following projection tensor V (bond dimension \ = 2) is 



the true ground state of the system: V 



TT 



1,2 



1, and 



T^i = = — 1 (the negative sign originates from the mi- 
nus sign between spin up and down components in Eq. ([3]l); 
similarly, for J\ > 0, the limit case locates at region II, and 
the ground state on odd bond is 



\<t> af (2i ~ 1, 2i)) = i| Ui-xhi - Ui-ihi), 



(4) 



with energy —J\ — J 2 , and the tensor T in ground-state MPS is 

-1. In addition, for 
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both cases, A=-^/isa2x2 diagonal matrix with doubly 
degenerate values. 

In the practical iTEBD projection process, the MPS wave- 
function is usually organized as two-period, i.e., it consists of 
two types of tensors T a and Tf,, 



N' 



= Tr(Q A.T™-- 1 A b T™*%.., m 2i . 1 ,m 2i 



(5) 

where A a (Af,), is Xa(b) x Xa(b) diagonal matrix on the cor- 
responding bond. The exact MPS expressed with one tensor 
r can also be rewritten with T a and Tf,. Because the spins 
on even bonds are unentangled, bond dimension Xa — 1 an d 
A a = 1, the other bond dimension Xb = 2, with diagonal ma- 
trix Ah = 

are (T b j[ A 
while for J\ 

{Ta)^ = (T a )i, 



When Ji < 0, the nonzero tensor elements 

- 1, (T 6 )k = -1, and (T a )l tl = (T a )f i2 = 1; 
0, (TOt,! - 1, (T b )t A 
1. 



> 



1, and again 



11,1 — \* a J 1,2 

Besides regions I and II, there also exist exact MPS ground 
states in regions III and IV. In Fig. Q] those parameter points 
are along the line J 2 = 0. Owing to the absence of quantum 
fluctuations, the model reduces to a classical Ising model with 
alternating couplings J\ and L\, and the ground state is direct 
product state, i.e., MPS with bond dimension \ = 1. When 
Ji > 0, one ground state spin configuration is illustrated in 



Fig. |2](b), and the MPS is two-period, with (T a ] 



1,1 



1 and 
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(T b )\ 1 = 1. While for J x < 0, one spin configuration is 
shown in Fig. |2](a), and the MPS is four-period, i.e., consists 
of four T tensors. Ignoring the bond indices owing to % = 1, 
the nonzero elements are T^ n _ 3 = T^ n _ 2 = 1, and T^ n _ 1 = 
T^ n = 1, where An — i is used to mark the lattice site, and 
n = 1,2, N'/2 (N' is assumed to be even number). 

Apart from the above limiting points, the MPS ground state 
can not be written down generally, however, we can adopt the 
ordinary perturbation theory to argue that MPS description is 
still a very nice ground state approximation. Firstly we take 
regions I and II as examples, and consider the lowest excited 
odd bonds as following, 32 

1^/(2* -1,2*)) = -^=| tu-iU+hi-M, (6) 

and 

KM2i- 1,2*)) = -j=\ t2i-i4-w +4*i-it2i), (7) 

with bond energy J\ + J2 and — J\ + J2, respectively. Where- 
after, we denote the one-particle excited state 

|S(i)) = |...^/(2i-l,2*)^/...), (8) 

for J 1 < 0; and 

\E(i)) = \...<t> af Tp af (2i - 1, 2t)<f> af ...), (9) 

for Jx > 0. That is, one odd bond (2i — 1,2*) is in state 
\il>f( a f))> while others remain in \4>f( a f))- It i s straightfor- 
ward to verify that the transition matrix element of pertur- 
bation operator between \E(i)) and zeroth-order exact MPS 
ground state vanishes, i.e., (^MPs|£i^fjCrf ' +1 |i£(*)) = 0- 
This fact suggests that lowest one-particle excited states do 
not affect the MPS wavefunction in the first-order approxima- 
tion, and only the multi-particle excited states or higher order 
perturbations will modify it. Further analysis shows that the 
perturbation term Lx cr | 7 ' <T 27+i w iH move the "particle" along 
the chain, i.e., (E^lLia^i&ii+i + 1)) = -^1> so tne one_ 
particle excitation dispersion up to the first-order approxima- 
tion can determined as^ 

4'" = 2J 2 + 2L lC os(fc), (10) 

in which k = -n + - jfi,Tr, and \k) = 

^= YaLi e lkl \E(l)). This dispersion suggests that excita- 
tion gap of the system is nonzero in the phase diagram except 
for the line J2/L1 = 1 (gapless at k — it). 

For regions III and IV, where the term ^ffi-i 17 ^ i s re ~ 
garded as perturbation, the same conclusion can be drawn af- 
ter similar arguments, i.e., single-particle excited states will 
not modify the MPS ground state up to the first-order single- 
particle perturbation, and MPS is also a very nice approxima- 
tion in these two gapped regions. Which is different, in this 
case the excited particle is revealed to be a moving "domain 
wall" instead of a single excited odd bond, and the dispersion 
relation can be verified as e l ^ l ' 1Y — 2L\ + 2 J2 cos(fc). 



C. iTEBD and imaginary time projections 

Beyond the perturbation arguments, imaginary time projec- 
tion technique iTEBD is employed to accurately determine 
the variational MPS wavefunction. 27 To be concrete, the vari- 
ational ground state \^ g ) (in MPS form) can be obtained by 
acting the imaginary time evolution operator exp(-f3H) on an 
arbitrary initial state |^o)- The operator exp(—(3H) is ex- 
panded through Suzuki-Trotter decomposition as a sequence 
of two-site gates JJ^ l ' l+1 ^ = exp(— r/i^.i+i), where hi.i+i is 
the local bond Hamiltonian, and r means small Trotter step 
length. In the limit j3 — > 00, the resulting wave function 
exp(-/?_ff)|^E'o) will converge (or be very close) to the ground 
state of H. Fig. [3] illustrates calculation errors compared with 
the exact solutions. Some typical parameters including criti- 
cal and noncritical points are concerned. The errors converge 
rapidly with enhancing x, in noncritical regions very accu- 
rate results can be obtained even with the smallest nontrivial 
bond dimension \ — 2,4, which convince us that MPS de- 
scription of the present system is not only adequate but also 
highly efficient and accurate. In practical implementations, 
the convergence of results with different bond dimension \ 
has always been checked, and for most cases up to \ = 40 
is quite enough. The total number of iterations taken is about 
10 5 — 10 6 . We first start with a step r = 10~\ and then di- 
minish it to t = 10~ 8 gradually. Whenever t is small enough, 
this procedure would bring the MPS to its canonical form, 
which would be useful for calculating the entanglement en- 
tropies, as well as the local observables including energy per 
site and local magnetizations, etc. 

During the iTEBD process for two-period MPS (Eq. (0), 
only four tensors (T a , Tb, A Q , and Ah) are involved and up- 
dated in each iteration step. In order to capture more sym- 
metry broken phases with larger unit cell, sometimes we need 
four-period MPS which includes eight different tensors (T a , 
Tfc, T c , Td, A a , Ah, A c , and A^). For example, region III in 
Fig. [TJis verified as a stripe AF ordered phase (one such or- 
dered spin configuration is illustrated in Fig. 12(a)), and this 
stripe AF order can be well described with four-period MPSs, 
but not two-period ones. 

D. Quantum Entanglement and Fidelity 

Quantum entanglement has close relationship with QPTs 
in many-body systems, 34 and much effort has been devoted to 
studying the quantitative description of entanglement in quan- 
tum systems^di In order to describe the QPTs in the EQCM, 
the von Neumann entropy SVn is adopted as a bipartite entan- 
glement measure. 40 When the MPS is gauged to its canonical 
form, i.e., we can cut an arbitrary bond in the system, and 
obtain a Schmidt decomposition as, 

|*) = f>^)A Q |<I>f). (11) 

a = l 

Here, |<£>^) (1$^)) represent the orthonormal bases of subsys- 
tem to the left (right) of the broken bond, and A is a diagonal 
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matrix. Correspondingly, the canonical MPS would satisfy 
the following two equations, 



m a 

EE( r *)™^ A ' r ™"^ = 

m ft 



(12) 



The superscript * means complex conjugate. It is easy to 
check that in the above limit cases L\ = 0, the exact MPSs 
satisfy Eq. (fT2l . and they are thus in canonical form. Given 
the MPS in its canonical form, bipartite entanglement of half 
chain (S'haif) can be directly read from diagonal matrix A (see 
Eq. CD), 



'half 



-Tr(A 2 log 2 A 2 ) 
-f»g 2 A 2 . 



(13) 



Notice, for two-period MPS, we can define two different bi- 
partite entanglement measures S 2 ;_i. 2 i = — Tr(A 2 log 2 A 2 ) 
and S2i,2i+i = — Tr(A 2 log 2 A 2 ), on odd and even bonds, 
respectively. Besides S'haif, people are also interested in the 
block entanglement, which is defined as follows, 

Sl = -Tr[p L log 2 (/? L )] = -Tr[p onVL log 2 (ponv L )], (14) 

where is the reduced density matrix of the L spin sys- 
tem, and envL means the environment (rest of the chain). Sl 
characterizes the entanglement between L adjacent spins and 
the environment. In practical calculations, the density ma- 
trix Penv L is supported by \ Schmidt bases, and employed 
to calculate the block entanglement entropy for arbitrary spin 
portion L. We would like to stress that the two kinds of entan- 
glement measures are both von Neumann entropies just that 
the bipartition happens to be between the left and right halves 
in the first case and between a block and the rest in the second 
case. 

Except for the entanglement measures, fidelity per site / is 
also used to detect the QPTs, 42 which is defined as 



/ = lim 



rcf/ 



N- 



N 



(15) 



|\&) is the ground-state wavefunction of the present system, 
and l^rcf) is a reference state, / indicates how fast the overlap 
of two distinct state decays to zero with increasing the length 
of the chain. The bifurcation and singular points of / can be 
utilized to locate the QPTs42i2 

Remarkably, the von Neumann entropy and the fidelity per 
site / can be conveniently obtained in the framework of MPS. 
Therefore, we will adopt them, along with the energy, magne- 
tization, and nearest-neighbor correlators, to study the phase 
transitions of ID EQCM in the following sections. 



m. QPTs in the one-dimensional EQCM 

E. Ground State Energy, Entanglement Entropy, and Local 
Order Parameter 

Firstly, we consider the QPTs along the line J 2 = 2.0 x 
(1.0 — Ji) (dash line (a) in Fig. [TJ- As the phase diagram 
illustrates, with increasing J\, the system should undergo two 
sequential QPTs: one first-order QPT from region I to region 
II and then the second-order one from region II to region IV. 
The bipartite entanglement entropies S2»-i,2« an d S 2 i j2 i+i are 
plotted in Fig. H](a). From Fig. |H(a), it is clearly seen that 
there exists only one singular point J\ = 0.5 (and J 2 = 1.0) 
where a second-order QPT takes place. From Fig. |4](b), an 
energy level crossing happens at J\ = 0, which indicates that 
a first-order QPT should occur there. However, as shown in 
Fig. |4] (a), the bipartite entanglement changes continuously 
across the first-order QPT. 44 Therefore, the first-order QPT 
in EQCM is missed by the entanglement measurement Shaif • 
Notice that the adiabatic continuations are plotted with dashed 
lines in Fig. |4](b), which illustrate the adiabatically evolved 
states from the left (or right) of the transition point, explicitly 
revealing the nature of level crossing at the first-order QPT. 13 

Next, we pay attention to the ground-state energy on odd 
and even bonds (denoted as e 2 i-i i2 ; and enpi+i, respec- 
tively) and their first-order derivatives (see Fig. 0(a) and (b)). 
We find that, the first-order QPT at J\ = can be detected 
by the energy level crossing of the odd bond energy (Fig. [5] 
(a)) or the discontinuous behavior of its first-order derivative 
(Fig. |5](b)). Furthermore, the singular behavior of the first- 
order derivatives (of both e 2 i_i j2 i and e 2 i. 2 i+i) at J\ = 0.5 
indicates the occurrence of the second-order QPT. According 
to the Feynman-Hellmann theorem 



de dH(X) 



(16) 



where A is a tunable parameter in the Hamiltonian. One 
can speculate that the first-order derivative of bond energy 
is in fact a second-order derivative of site energy e. Take 
even bond energy e 2 i. 2 i+i as an example, de2i,2i+i / dJ\ — 
d 2 e/(dJidLi), and it is thus expected to show singular be- 
haviors around the second-order QPTs (as Fig. 0(b) shows). 

On the other hand, de^i-ifiil dJ\ = (c|i-i (T 2i) + 
J2d(a-2 i _ 1 (J2 i ) /dJi, and the short-range correlators 

(&2i-i cr 2i) an< ^ ( cr fi-i <7 2i) on °dd bonds are calculated and 
shown in Figs. 0[7] From Fig. we find that the short-range 
correlation (&2i-i a 2i) I s +1 m region I, but abruptly changes 
into —1 as entering into regions II and IV. So, the first-order 
QPT takes place with a sign change of two-site correlation 
function (cf^crfi) on odd bonds, and causes a discontinuity 
in de-2.i-\,nl dJ\ curve. In Fig. |7]( a ), the two-site correlator 
( <7 2i-i f7 2i) on °dd bonds behaves continuously, and the 
divergent peak of its first-order derivative (Fig. Q (b)) can 
signal the critical point ( J\ = 0.5). Consequently, derivative 
de^i-ifli/dJx also has a divergent peak at J\ = 0.5. At last, 
it is worth noticing that, although derivatives de^i-ifii/ 'dJ\ 
and de2i t 2i+i/dJ% are both divergent at J\ = 0.5, they have 
different signs and cancel with each other, the first-order 
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FIG. 4: (Color online) (a) Half-chain entanglement entropy on odd 
bond S-2i-i,2i and even bond 5*21,21+1 . (b) Ground state energy per 
site, two dotted lines represent energies of the adiabatically evolved 
states from left and right sides of the first-order QPT point. 




J/L, J/L, 

FIG. 5: (Color online) (a) The odd (even) bond energy e& and (b) 
their first-order derivatives along the line J2 = 2(1 — Ji). 



derivative de/dJi — d(e2i-i,2i + e2i,2i+i) / dJi is continuous 
across the second-order QPT point. 

To attain comprehensive understandings of the QPTs in 
EQCM, we also compute the local magnetizations \(a x ) \ and 
I (a z ) I (shown in Fig. [8]), whose values are independent of the 
sites owing to the translational invariant MPS. It is observed 
that across the critical point J\ = 0.5, the EQCM goes into 
a region (IV) with nonzero |(er 2 )|, and their values have dif- 
ferent signs on odd and even sites, i.e., staggered magnetiza- 
tion, and thus region IV can be recognized as a semi-classical 
Neel phase. Thus, the magnetization |(cr z )| (to be more strict, 
staggered magnetization M^ eel = §|{cr|i— 1 — °2i)D can be 
recognized as the local order parameter characterizing region 
IV in Fig. Q] however, the regions I and II are both disor- 
dered phases, and can not be distinguished by the local order 
parameters. 45 On the other hand, Fig. [8] reveals that | (cr^) | 
vanishes in regions I, II, and IV on either odd or even sites. 

In order to discuss the bipartite entanglement behavior 
across the multicritical point ( J\ = 0, J2 = 1.0), we consider 
the line J2 = 1.0 — J\. Along this line, with increasing J\, 
the ground state of EQCM will go from region I into region 
IV through the multicritical point. The odd and even bond bi- 
partite entanglement measures are plotted in Fig. [9] (a). It is 



FIG. 6: Nearest neighbor correlator {crfi-i "!*)- 
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FIG. 7: Short-range correlation (o"2i-iO"2i) ( a X and its first-order 
derivative (b). 



noticed that the QPT can be recognized by the sharp peaks of 
S2i-i,2i and S2i,2i+i, which confirms it as a quantum critical 
point. It should also be mentioned that, except for the sim- 
ilar second-order QPT character in entanglement measure, a 
distinctive ground-state energy level crossing is also clearly 
shown in Fig. 0(b), where adiabatic continuations are again 
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FIG. 8: (Color online) Local magnetizations (cr^ ) | and |(<r*)|. 
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FIG. 9: (Color online) (a) Bipartite entanglement measures Ssi-i.ai 
and £21,21+1 along the line J2 = 1.0 — Ji. (b) Ground state energy 
per site e, the dotted lines represent energy of adiabatically continued 
states. 



FIG. 10: (Color online) (a) Bipartite entanglement on odd bond 
£21-1,22 and even bond £21,2^+1 along line J2 = 0.8(1.0 — Ji). 46 
(b) Ground state energy per site e (the dotted lines are adiabatic con- 
tinuations). 



employed to verify this conclusion. Therefore, both the first- 
and second-order QPT features at this multicritical point are 
revealed by our calculations. Besides, the local magnetization 
\{<J X ) \ and \(<t z ) \ in regions I and IV are also evaluated (not 
shown for the sake of space), and similar behaviors from dis- 
ordered region I to Neel ordered region IV as discussed above 
are again observed. 

Lastly, we consider the QPTs along the line Ji = 0.8 x 
(1.0 — Ji). With increasing J±, the ground state of EQCM 
will go from region I to region III, and then enter into region 
IV. The bipartite entanglement on odd and even bonds, and the 
ground-state energy per site e are plotted in Fig. [I0](a) and (b), 
respectively. We find that, although the second-order QPT at 
Ji = —0.25, Ji = 1 is signaled by a singular peak of the 
entanglement entropy, the first-order QPT at Ji = with dis- 
tinct ground-state energy level crossing (Fig. [I0](b)) is again 
missed by the entanglement measure (Fig. H0](a)). However, 
as shown in Fig. QT|(a) an d (b)> the bond energy and their first- 
order derivatives are able to capture all the QPTs. In addition, 
magnetization is calculated and shown in Fig. [12] nonzero 
| (a z ) | is found in region III and IV, and | (a x ) | vanishes along 
the whole line on either odd or even sites. 

Moreover, in Fig. Q~2] although the magnitude of |(er 2 )| 
change smoothly through the phase transition point Ji = 
0, Ji = 0.8, the magnetic order is quite different in region 
III with that of region IV. Calculations indicate that, the cor- 
relators (cfi-i^fi) = 1 in region III, show distinct difference 
with those in Neel phase (region IV), where (cfi-i "!-*) = — 1- 
In fact, the magnetic order in region III is four-period stripe 
AF order, quite different from Neel order in region IV. In Neel 
phase, the spins are arrayed in "up-down-up-down" pattern 
(see Fig. 12(b)); while in stripe AF phase, they are in " up-up- 
down-down " arrangements, one typical spin configuration of 
such phase is illustrated in Fig. |2](a). The stripe AF order in 
ID EQCM was previously proposed in Ref. 23 with finite- 
size calculations by Lanczos method, and it is confirmed here 
by our results directly in the thermodynamic limit. 

As the Neel order parameter M^ ccl defined in region IV, 
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FIG. 11: (Color online) (a) Ground state energy of odd and even 
bonds and (b) their first derivatives. 



Stripe = ~ ^4n-i)l- where n = l,2,...,/V'/2, 

can be defined as the local order parameter in region III. Then, 
Mj^ ool is nonzero in region IV, and vanishes abruptly in region 
III; while the reverse is Af s z tripc , which appears in region III, 
and drops to zero in region IV. Therefore, the first-order QPT 
between regions III and IV can be recognized by evaluating 
local order parameters, quite different with the transition be- 
tween regions I and II discussed above. 



F. Block Entanglement Entropy 

Besides the half chain entanglement, the block entangle- 
ment entropy Sl are also calculated, which provides a mea- 
surement of the amount of entanglement between L adjacent 
spins and the rest of the system (environment). With MPS 
wavefunction, we are able to obtain the reduced density matrix 
of the environment supported by the bond bases, and hence 
can calculate the Si, with length L up to several hundreds 
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FIG. 12: (Color online) Local magnetization | {cr 31 ) | and |(cr z )| (J2 = 
0.8(1 — Ji )). The magnetic order is stripe AF in region III, and Neel 
AF in region IV. 
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FIG. 13: (Color online) (a) Entanglement entropy Sl (with block 
size L — 4) along line J2 = 1.2(1.0 — Ji) and (b) Sl saturates 
rapidly with increasing L in the vicinity of the first-order QPT line. 



of sites at ease. The block entanglement entropy (Sl) with 
L — 4 along the line J 2 = 1.2 x (1.0- Ji) is plotted in Fig. [13] 
(a). With increasing J±, two sequential QPTs will take place: 
one first-order QPT from region I to region II and the other 
second-order QPT from region II to region IV. However, from 
Fig. [I3](a), we find that only the second-order QPT at J 2 = 1 
can be detected by the peak of the block entanglement entropy 
Sl, the first-order QPT (at Ji = 0, J 2 = 1.2) between phase I 
and phase II is missed again. Sl continuously approaches the 
same value whether J\ — > 0~ or J\ — > + with fixed J 2 . In 
Fig. [T3](b), for the non-critical ground state, when block size 
L increases, the block entropy Sl enhances and quickly be- 
comes saturated, well satisfying the entanglement area law. 30 
These observations on the block entanglement are consistent 
with those proposed in Ref. |22[ Therefore, the entanglement 
measures, including block entanglement entropy Sl and half 
chain entanglement Shaif, are indeed not able to detect the 
first-order QPTs in ID EQCM. 



Next, the scaling behavior of the block entropy Sl on the 
second-order QPT line J 2 = 1 is investigated. As shown in 
Fig. El the block entropy Sl exhibits divergent behavior with 
increasing block size. As derived in Ref. |47l in a 1 + 1 di- 
mensional conformal field theory, the entropy of a subregion 
of length L reads 

Sl = ^log 2 (L) + k, (17) 

with a coefficient given by the holomorphic and antiholomor- 
phic central charges c and c of the theory. From Fig. [14] we 
find that the divergent Sl on the second-order QPT line can be 
well fitted by Sl = ^log 2 (L) + 0.5202, with central charges c 
= c=l/2, i.e., the Sl displays a logarithmic divergence on the 
second-order QPT line. Therefore, we disclose that the criti- 
cal behavior of EQCM can be described by a free fermionic 
field theory, 34 with central charges c/ = Cf =1/2. 



G. Fidelity Calculations 

Except for the entanglement, the fidelity measure defined 
in Eq. ([151) is also utilized to study the QPTs in EQCM. 
Facilitated with MPS framework, it is straightforward that / 
can be obtained by evaluating the maximum eigenvalue of the 
transfer-matrix defined as 

m 

in which T (along with A) represents the reference state. Con- 
sidering the multi -period MPS wavefunctions (period 2 for re- 
gions I, II, and IV, and period 4 for region III), we slightly 
modify it and define the fidelity per unit cell, which is the 
maximum eigenvalue of the transfer-matrix defined in a unit 
cell. For instance, the transfer matrix of two-period MPS can 
be defined as following, 

^t,y 7 = E (A«)«'(r;)^7 1 (A i V(7?)^ 

/3',/3,m2i_i,m2» 

x (A a ) a (T a )^(A b ) p (T b )^, (19) 

which is a Xa x xl matrix. The transfer matrix of four-period 
MPS can be similarly written down. 

In Fig. [131 the results of fidelity per unit cell are present 
(the MPSs are generally set as period 4) along three different 
lines, J 2 = 2(1 - Ji), J 2 = 1 - Ji, and J 2 = (1 - Ji)/2, 
respectively. In Fig. Q3J (a), the line traverses regions I, II, 
and IV, and in Fig. Q3J (b), regions I and IV are involved. 
During these calculations, the ground state of Hamiltonian 
Eq. (Q]i with parameter J\ = 1, J 2 = is set as the ref- 
erence state (i.e., an Ising AF state). Owing to the sponta- 
neous Z2 symmetry breaking in the Neel phase, / shows bi- 
furcation behaviors in Fig. [I5j(a) and (b), and the bifurcation 
points locate the second-order QPTs. Besides, the first-order 
QPTs can also be recognized from the discontinuities in fi- 
delity curves. It is worth noticing that the results in Fig. [15J 
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FIG. 14: (Color online) Scaling of the block entanglement entropy 
Sl on the second-order QPT line J2 = 1.0, the solid line is fit to 
numerical data. 



(b) again reveal multicritical properties of the transition oc- 
curred at Ji = 0, J-2 = 1, i.e., the discontinuity of / indi- 
cates a first-order QPT, while the bifurcation phenomenon re- 
veals second-order QPT character. In Fig. [15] (c), we choose 
Ji = —0.5, J2 = 0.75 as the reference point, the bifurcation 
at Ji = —1 indicates second-order QPT between regions I 
and III, and the discontinuity at J\ — suggests the first-order 
QPT between stripe AF and Neel phases (regions III and IV, 
respectively). 



H. Entanglement Spectrum, Dual Transformation, and String 
Order Parameters 



Through previous analysis in subsection PEl it is uncovered 
that no local order parameter can be utilized to distinguish the 
two disordered phases regions I and II in Fig. Q] as well as to 
detect the first-order QPTs between them. In this subsection, 
the non-local string order parameters in regions I and II are 
computed and discussed. 

In Fig. [16] several typical entanglement spectra of ID 
EQCM are shown, which exhibit the eigenvalues of the re- 
duced density matrix of half-infinite chain by dividing the 
system via any bond. For the canonical MPS, entanglement 
spectrum can be recognized as the diagonal elements of A 2 in 
Eq. ( fTTb - For the present two-period system, it is free to cut an 
even or odd bond, thus we have two entanglement spectrums 
(A 2 and A 2 ) for a single parameter point. In Fig. [16] (a), (b), 
and (d), non-critical points are concerned, and the eigenvalues 
decay roughly exponentially; while in Fig. [16] (c), for critical 
point, the entanglement spectrum decays much slower, and in 
some algebraical way. Another distinct feature is the doubly 
degenerate A 2 for disordered phases I and II (Fig. [16] (a) and 
(d)), which implies the existence of the non-local string order 

AO 

parameters. 
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FIG. 15: (Color online) Fidelity per unit cell along three differ- 
ent lines are present in (a), (b), and (c). The discontinuities of the 
curves indicate the occurrence of first-order QPTs, and the bifurca- 
tion phenomena manifest the spontaneous Z2 symmetry breaking, 
where second-order QPTs take place. Four-period MPS is adopted 
during the calculations. Notice that the reference state in (c) is dif- 
ferent from that in (a) and (b), see the text for more information. 



exist two topological distinct disordered phases for JijL\ > 
1 and J2/L1 < 1, and the phase transition between them 
(at J2/L1 = 1) is disclosed as a topological QPT, 49 char- 
acterized by non-local string order parameters. Other than 
this disordered line, it is an interesting question that whether 
the string order parameters in regions I and II still exist or 
not. To accomplish this task, standard Kramers-Wannier dual 
transformation 50 is employed to map the present model to the 
quantum-Ising system. The dual mapping of each terms in 
Hamiltonian Eq. are as followings (here we adopt the for- 
malism introduced in Ref . 5 1 , and a permutation of even and 
odd bonds is taken before dual transformation), 



- / 2<7 2l C7 2i+1 -> -J 2 CT i (T i+1 , 



and thus the Hamiltonian is as 



iV' 



IT; T. 



i+1 



(20) 



(21) 



where a and t are Pauli matrices on the dual lattice. Dual 
Hamiltonian Eq. d2Tb can be regarded as two decoupled Ising 
spin chains (couplings —J\ and — J2, respectively), 52 and the 
a chain is under transverse field {L{af term). There may exist 
two types of spontaneous long-range orders, i.e., (o^er 2 ) and 
Previous studies suggested that along the line J\ — there (t^t 2 ), which can be mapped back to the original system as 
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the following non-local string order parameters, 

(- 1 Y l ~ k ( T k T n) -> ( cr 2fe cr 2fc+l-" cr 2n-2 cr 2«-l)> 

{-l) n - k {alal) {al k al k+1 ...al n _ 2 al l _ 1 ). (22) 

The two types of a operator strings can be denoted as O zz (ri- 
le) and O xx (n — k), respectively. Owing to the absence 
of transverse field on t spins in the dual model Eq. OTT) . 
O zz (n — k) is always nonzero in the whole phase diagram. 
To be specific, it is found that O zz (n — k) = 1 for region I 
(and also region III), while O zz (n — k) = (— l) ,l ~ fe for region 
II (and IV). This is owing to that in the dual model there exists 
ferromagnetic long range order ((t^t*) = 1) for J\ > 0, and 
AF long range order ((t£t*) = (— l) n ~ fe ) for Jj < 0. On the 
other hand, this conclusion can be easily verified by noticing 
the nearest neighbor correlators (fffi— l^Ii) = •"■ f° r re gi° n 1 
(III) and —1 for region II (IV), which can also be regarded as 
good quantum numbers for ground states. 

The behavior of the other string order parameter O xx (L) 
(L = 2(n — k) is the number of sites in the string) is more in- 
triguing, and the numerical results are shown in Fig. [17] The 
inset of Fig. [l7]shows that the O xx (L) converges very rapidly 
with L (except for points in the vicinity of second-order QPT 
line J2 = 1). The converged O xx (oo) monotonously de- 
creases with enhancing the parameter J\ (and hence decreas- 
ing J 2), and changes continuously through the first-order QPT 
line Ji = 0, vanishing immediately after crossing the second- 
order QPT line J% = 1. The asymptotic behavior of O xx 
near line J 2 = 1 can be predicted by the dual spin correlation 
function,^^ as |O xa (oo)| ~ (1 — Jf) 1 / 4 , which can be well 
verified from the fitting in Fig. [TTl 

Therefore, the above investigations uncover that in regions 
I and II, the string order parameters O xx and O zz are nonzero, 
which reveals the hidden Z2 X Z2 symmetry breaking in the 
EQCM system. In addition, O zz can be used to distinguish 
two disordered phases and detect the first-order QPTs between 
them, while O xx changes continuously through the transition 
line Ji = 0, and vanishes at critical line J 2 — 1. On the 
other hand, it is reported in Ref. 54 that the nonzero string or- 
der O xx in disordered region is robust even under some finite 
external magnetic fields h < h c (below the critical field). 
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FIG. 16: (Color online) The entanglement spectra of several noncrit- 
ical (a, b, and d) and critical (c) points. In (a) and (d), A5 's are doubly 
degenerate. 
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FIG. 17: (Color online) The nonlocal order parameter O xx , inset 
shows the behaviors of O xx (L) with portion length L (shown for 
every two sites). The converged value O xx (oo) is nonzero in dis- 
ordered regions I and II, and vanishes in ordered phases III and 
IV. The fitting lines illustrate the asymptotic behavior \O xx (oo)\ ~ 
(1 — Jf ) in the vicinity of the second-order QPT. 



IV. Specific Heat Curves 

Besides the ground-state properties, in this section the 
LTRG method 55 is employed to investigate the finite tem- 
perature properties of ID EQCM. LTRG adopts the iTEBD 
technique for contracting the transfer-matrix tensor network, 
and can accurately (and efficiently) obtain the thermodynamic 
quantities including free energy, energy, susceptibility, and 
specific heat, etc. In Ref. |55_l LTRG method has been applied 
to calculate the isotropic XY model and achieved very accu- 
rate results. In order to verify the applicability and accuracy of 
LTRG for the anisotropic cases (for the present EQCM, there 
exist strong anisotropics in spin couplings), the specific heat 



curves of anisotropic XY model with Hamiltonian 

<i,j> 

are calculated, and shown in Fig. [18] The results of LTRG 
show perfect coincidence with the exact solutions in Ref. H^. 

In Figs. [19] [20] the specific heat (C) curves of ID 
EQCM are present, and close attention is paid to their low- 
temperature behaviors, which reveal the low-energy excita- 
tion features of the system. In Fig. [19] (a), the specific heat 
curves are evaluated along the critical line J2/L1 = 1. It is 
observed that with gradually decreasing | J\ |, there appear low 
temperature sub-peaks moving towards T = 0, and disappear 
when Ji = 0. In Fig. [19] (b), the low temperature T parts 
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FIG. 18: (Color online) Specific heat curves for anisotropic XY 
model. The solid and dashed lines represent the exact solutions, and 
the scatters are the LTRG results which show perfect agreements with 
the lines. The couplings J x and J y are defined in Eq. J23b . Jy = 1 
is set as energy scale here. 




FIG. 19: (Color online) (a) Specific heat (C) curves along the critical 
line J2 /Li — 1. (b) The low temperature sections of C. The specific 
heat curves with the same absolute values | Ji | (but different signs) 
almost coincide with each other at low temperatures. The exact solu- 
tions (see Ref. 57) are also plotted with lines, with which the LTRG 
results show very nice agreements. 



of the C curves are magnified, and the linear relations with T 
are clearly shown, which can be ascribed to the gapless low- 
energy excitations along the critical line of EQCM. 

Along the parameter line JijL\ = 1.5 , the C curves (ver- 
sus temperature) are illustrated in Fig. [20] where sub-peaks 
also appear and similar movement behaviors are again ob- 
served in Fig. [20] (a). It is worth noticing that there exist 
excitation gaps along the line J2/X1 = 1.5, as revealed in 
Fig. [20] (b), where the low temperature C with Ji = is 
shown to decay exponentially. Furthermore, from Fig. 1201(b) 
(judging from the slop of the C curves in the log-log plot), 
it is found that the excitation gap A(Ji) tends to zero when 
the parameters approach J\ = line from both sides, but the 
Ji = 0, J2/L1 = 1.5 point itself is far from gapless. 




FIG. 20: (Color online), (a) Specific heat (C) curves along the line 
J2/X1 = 1.5. (b) The low temperature parts of C, log-log plot 
reveals the exponential decay explicitly, and the slope of the (nearly) 
straight part of C is intimately related with the excitation gap. In 
(b), the lines of Ji = ±0.5 (as well as Ji = ±0.8) almost coincide 
with each other, and Ji = ±0.8 lines are also very close to the line 
Ji = 0. 

V. Summary and Outlook 

I. Subsequent Problems 

The exact soluble ID EQCM provides an ideal playground 
for performing calculations and testing MPS-based algo- 
rithms, remarkable accuracy and perfect accordance with pre- 
vious analytical results have been achieved. Nevertheless, we 
would like to stress that the power of MPS-based numeri- 
cal methods is their accessibility to more complex problems 
which do not permite exact solutions. Among others, we re- 
gard the following two extended compass models particularly 
interesting, 

N' 

HxYZl = ^ Jxx^n-l^n + J yy a 2n~\ a \n + ■4z cr 2n' 7 2 n+1 1 

71=1 

(24) 

and 

jv" 

#XYZ2 = ^2 J ^ a lri-2 a tn-\+ J Vy a ln-\ a ln + Jzz a ln a ln+\i 
ra=l 

(25) 

where N' = N/2, and N" = N/3, N is the total site number. 
The first one is two-period, with J xx and J yy couplings on 
odd bonds and J zz couplings on even ones, the second model 
is three-period, with J xx , J yy , and J zz on three different types 
of bonds, respectively. These two models are more complex 
than EQCM in Eq. (fl}, while are still very basic ones. Owing 
to the existence of three non-commuting spin coupling com- 
ponents in the Hamiltonian, they can not be diagonalized by 
simply taking fermionic transformations, and are expected to 
show more interesting QPTs in their phase diagrams. In the 
first model Hxyzi, compared with Hamiltonian Eq. (Q]), the 
J zz couplings on odd bonds are replaced with J yy couplings, 
so the expectation values of parity operators cfi-i '!?; on 
bonds are no longer good quantum numbers. Some prelimi- 



12 



nary results are obtained by iTEBD calculations, which reveal 
that there also exist first- and second-order QPTs, as well as 
multi-critical points in the phase diagram of model -f/xYZi- 
A distinct difference between the phase diagram of -Hxyzi 
and Fig. Q]is that the Neel and stripe AF zones in the present 
EQCM are extending along the J\ axis to the infinity, while 
for the former case Hxyzi, they are confined in a finite re- 
gion. More details about the ground-state phase diagrams and 
QPTs in these two EQCMs will appear elsewhere. 

J. Conclusions 

Employing MPS wavefunction, and with the aid of the 
related algorithms iTEBD and LTRG, we investigated the 
ground-state properties and QPTs, as well as specific heat 
curves, in the ID EQCM. 

Our calculations, including energy per site, bond energy, 
entanglement entropy and local magnetizations, validate the 
phase diagram proposed by previous works. Four different 
phases are identified in Fig. [Tj including two disordered re- 
gions I and II, Neel ordered phase (region IV), and a stripe 
AF phase in region III. 

The second-order QPTs along J2/L1 — 1 line can be de- 
tected by the singularities of entanglement entropy, as well as 
the derivatives of bond energy. The first-order QPTs along 
Ji = are however indeed missed by entanglement measures 
according to our calculations. Furthermore, at the multicrit- 
ical point ( J\ = 0, JijL\ — 1.0), besides the second-order 
QPT feature revealed by entanglement entropy, a distinctive 
ground-state energy level crossing (observed by taking adi- 
abatic continuations) occurs. Therefore, at the multicritical 
point, there coexist both the first-order and the second-order 
QPT characters. Furthermore, a logarithmic divergent be- 
havior of block entanglement Sl on the second-order QPT 
line J2/L1 — 1 are observed, from which the central charge 



c = 1/2 is determined. 

Fidelity per unit cell is also used to investigate the QPTs, 
and it is disclosed that both the first- and second-order QPTs 
in the EQCM can be detected by identifying the discontinuous 
and bifurcation points in calculated fidelity curves. 

Moreover, the disordered regions I and II are found to pos- 
sess doubly degenerate entanglement spectra, as well as two 
types of nonzero string order parameters O xx and O zz . By 
taking dual transformations, it is revealed that the string order 
parameters reflect the hidden Z2 X Z2 symmetry breaking, 
and parameters O zz can be used to detect the first-order QPT 
between regions I and II. 

Subsequently, the specific heat curves have been studied via 
LTRG calculations, and low temperature linear behaviors are 
observed along the critical line J2/L1 = 1, while for J2 7^ 
L\, the exponential decay of C at low temperatures implies 
the existence of a nonzero excitation gap. 

In conclusion, the fidelity per unit cell is shown to be sen- 
sitive to detecting not only the first-order but also the second- 
order QPTs, while the entanglement measures can only de- 
tect the latter ones. In the phase diagram Fig. [T| there exist 
two symmetry broken phases in regions III (stripe AF) and 
IV (Neel) with different local order parameters, and two hid- 
den symmetry broken phases in regions I and II with nonzero 
string order parameters. 
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